arXiv:1505.01330vl [physics.acc-ph] 6 May 2015 


Periodic Poisson Solver for Particle Tracking 

M. Dohlus, Ch. Henning 

Deutsches Elektronen-Synchrotron (DESY), D-22607 Hamburg, Germany 

Abstract 

A method is described to solve the Poisson problem for a three dimensional source distribution 
that is periodic into one direction. Perpendicular to the direction of periodicity a free space (or 
open) boundary is realized. In beam physics, this approach allows to calculate the space charge 
field of a continualized charged particle distribution with periodic pattern. 

The method is based on a particle mesh approach with equidistant grid and fast convolution 
with a Green’s function. The periodic approach uses only one period of the source distribution, 
but a periodic extension of the Green’s function. 

The approach is numerically efficient and allows the investigation of periodic- and pseudo-periodic 
structures with period lengths that are small compared to the source dimensions, for instance 
of laser modulated beams or of the evolution of micro bunch structures. Applications for laser 
modulated beams are given. 
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I. INTRODUCTION 


The fundamental problem of the analysis of beam dynamics is to solve the equation of 
motion for a multi particle system (N ~ 10 3 • • • 10 9 ) in presence of external forces (caused by 
magnets, cavities etc.) and self forces. A particular difficulty is to determine the self-fields of 
the dynamic many-particle system. The coarsest approximation to solve this problem is the 
space charge or Poisson approach. It assumes that all particles are in uniform motion with 
identical velocity vector. Therefore the problem is equivalent to an electrostatic problem, and 
in principle one has to add for each observer particle the one over distance contribution from 
all other sources. This point-to-point interaction would allow the slightly better approach 
of individual uniform motion per particle, but its effort scales as N 2 . This is avoidable by 
a continualization so that the point particles are represented by a smooth density function. 
Additionally the continualization solves a second problem that is usually not in the focus 
of interest: the near scattering of particles that approach much closer than the averaged 
particle distance. A proper tracking with near scattering needs not only the point-to-point 
interaction but also a very fine time step related to the minimal distance between particles. 

Continuous density functions and fields are numerically represented either by superposi¬ 
tions of basis functions (finite element methods) or by mesh methods. In this report we use 
a three dimensional equidistant mesh for binning and compute the static potential by a fast 
convolution with a Green’s function as described in [1, 2], This Green’s function fulfills free 
space (or open) boundary conditions. Also the electric field is calculated on the mesh, by 
numeric differentiation. The field is interpolated to the particle positions and transformed 
into a collective force as it appears in the equation of motion. The numeric integration of 
the equation of motion is not described here. 

The Poisson-particlc-mesh approach is used widely and successfully in beam dynamics for 
instance in straight sections with constant energy, but also in particle sources or accelerating 
sections and even in sections with weak curvature. Sometimes the particle distribution has 


a fine sub-structure w 
plasma oscillations 


rich has to be taken into account, for instance the particles perform 
5j, a micro modulation builds up (micro-bunch-instability) or they 
interact in an undulator with a laser. Such high resolution computations need very fine 
meshes and high particle (or macro-particle) numbers so that they become numerically very 
expensive. 
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This numerical effort can be dramatically reduced when the approach assumes and utilizes 
periodicity: only the particles of one period have to be tracked and meshed. We describe 
such an approach where the periodicity of the source distribution is interchanged with the 
periodicity of the Green’s function. Except for the calculation of the Green’s function, the 
established method (Poisson-particlc-mesh with equidistant mesh and fast convolution) can 
be used completely unchanged. We consider periodicity into one direction. The direction is 
arbitrary and needs not to coincide with a mesh axis. Perpendicular to this direction free 
space (or open) boundary conditions are realized. 

The report is organized as follows. Section II is about general aspects of periodicity. 
In section III we describe the approach, first the standard method (Poisson, particle-mesh, 
convolution with Green’s function), then its modification with periodic Green’s function. 
In section IV we demonstrate the method with examples: the FLASH seeding section 6-81 
with different particle densities and modulations, and parasitic heating after the LCLS laser 
heater (i ll(l]. For comparison of field calculation methods, we did tracking simulations with 
Astra |nj and the non-periodic and periodic implementation of the method described in 
this report, called QField. 


II. PERIODICITY 

Some aspects of spatial periodicity and pseudo-periodicity are illustrated in Fig. Q] The 
particles in panel (a) are randomly distributed while panel (b) shows a periodically repeated 
set S p of random particles. Distribution (a) is caused by applying a periodic transport map 
to random particles that are uniformly distributed into the direction of periodicity. Strictly 
periodic distributions as in (b) are used to approximate pseudo-periodic distributions as in 
(a). 

For the strictly periodic approach we consider only a set S p of non-repeated particles. 
As it can be seen in panels (b-e) there are many sets S p that describe the same periodic 
distribution. It is easy to manipulate such sets by shifting particles by a multiple of the 
repetition vector r p . The sets in (c) and (d) are chosen so that all particles are in a volume 
slice, where the front plane corresponds to the back plane shifted by r p . These slices might 
have any orientation (or shift), but there are orientations so that the volume of a box (or 
sphere) that includes all particles is minimal. The distribution in (b) are particles from 
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one volume slice, that have been mapped by a periodic transport map. (Particles with 
nearly the same spatial coordinates, but with different momenta may be mapped to very 
different spatial coordinates.) For numerical reasons it might be advantageous to use sets 
without sudden change of density, as it can be seen in panel (e). This set has been composed 
similar as (c), but some particles close to the left boundary are shifted by r p (beyond the 
right boundary) and some particles close to the right boundary are shifted into the opposite 
direction. The probability of shifting depends on the distance to the boundary. 

To keep periodicity, the set S p has to represent a distribution that is periodic in six 
dimensional phase space X = [x, y, z,p x ,p y ,p~\. Periodically shifted particles are at X + nP , 
with P the vector of six dimensional periodicity. The required property for the transport 
T{X) —> Y is T(X + nP ) —> T(X) +nT(P). This is for instance fulfilled if only one spatial 
component of P is non-zero and T is a linear function of that component. 

Real distributions are always pseudo-periodic since they are finite. Sometimes only a 
small part of a (long) bunch behaves periodically, as for instance in Figs. [U O The phase 
space volume covered by the beam is usually small and the transport can be linearized, or 
at least locally (for one ore few periods) linearized. Therefore the Poisson approach can be 
used for distributions with finite extension in momenta space, although the uniform motion 
approach assumes identical momentum of all particles. In the same sense the periodical 
Poisson approach is applicable to distributions with non-zero momentum components of 
P , although this would mean particles with arbitrarily large momentum. An example is a 
bunch with energy chirp and micro modulation. 

The Coulomb force drops with 1 jr and the effect of a dipole perturbation drops with 
1/r 2 . Therefore space charge effects are very localized and non-resonant. Pseudo-periodic 
behavior can be observed even for distributions with only few periods, as it can be seen in 
the first example of section IV. 


III. METHOD 

A. Poisson Approach and Lorentz Force 

The electro-magnetic field caused by a set S of discrete point particles with charge 
position Ti(t) and velocity Vj(t) is estimated based on the assumption of collective uniform 
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motion (of all particles) with the velocity v c . The particle positions are represented by 
orthonormal coordinates r* = j_i + yi& j _2 + ^e c , wherein the 2 -axis points in the direction 
of motion e c = v c /||e c || and the directions e_|_i and e± 2 are perpendicular. After Lorentz 
transformation to the rest frame, the coordinates are r i c = q + y^e j _ 2 + 2 jj c e c with 
7 C = 1 /y/1 — f3% and /3 C = v c /c = ||v c ||/c. The electrostatic held may be computed from 
point-to-point interaction as 


Ej , 


1 

Atieq 


Y % 

jeSi 


Ti,c - ry c 



( 1 ) 


with excluded self interaction. Therefore Sj is the subset of all particles without particle 
i. The continualization replaces the particle distribution {qi, rj ]C } by a continuous charge 
density function p c (r) and the electrostatic held 


E c (r) = —V</> c (r) 


( 2 ) 


is the negative gradient of the electrostatic potential 


r = 


47T£n 


(IV X 


Pc(r') 


r — r 


Finally we transform back to the moving frame 


( 3 ) 


Ej — 7c (exi (exi • Ej iC ) + e± 2 (e ±2 • Ej iC )) + e c (e c • Ej iC ) (4) 

B, = c“ 2 v c x Ej (5) 


and calculate the Lorentz force 


Fj = e (Ej + Vj x Bj) = e (Ej + c 2 v, x v c x Ej) 


( 6 ) 


B. Particle-Mesh-Method 

A simple and robust continualization method is to approximate the discrete charges by a 
continuous distribution that is piecewise constant in the cells of an equidistant mesh. The 
mesh points r m (j,k,l) = jh x e _|_i + kh y e± 2 + lh z e c are in the middle of the cells and the 
charges per cell are q c {j,k,l). For simplicity we have chosen the direction of motion e c to 
coincide with one mesh axis. This has numerical advantages, as the resolution requirement 
in longitudinal direction is usually very different from that into transverse directions, but it 
is not necessary. The choice of the perpendicular directions ej_i, e± 2 is free. 
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The electrostatic potential caused by the charge density p = g c (0, 0, 0)/(h x h y h z ) in the 
cell which includes the origin, is 

0(r) = (7) 

4:7T£ 0 

with the Green’s function 


G(xe _li + ye ±2 + ze c ) = j 
This integral can be solved as 


f*h x /2 rhy/2 rh z /2 


dx'dy'dz' 


—h x /2 J-h y /2 J-h z /2 ■>/( X - X O 2 + (y- y') 2 + (* - z') 2 


( 8 ) 


fl: 


G(xe_Li + ye ±_2 + ^e c ) = ^ ijk H(x + ih x /2, y + jh y /2, z + kh z /2) , (9) 

with the anti-derivative 

TT , . x 2 yz y 2 zx z 2 xy 

H lx, y,z) =-arctan-arctan-arctan — 

2 xr 2 yr 2 zr 

+xyln(z + r) + yz\n(x + r) + zxln(y + r) , (10) 


and r = \Jx 2 + y 2 + z 2 . The asymptotic behavior 


G(r) G a ( r) = h x h y h z /r 


( 11 ) 


is obvious from Eq. (J8j). 

To compute the potential at all mesh points (j) m (j,k,l ), caused by all mesh charges, one 
has to perform the summation 


0m (j, k, l ) 


1 1 
47T£ 0 h x hyh z 


9c(/, k', l')g{j - j', k - 

j' ,k' ,1' 


k’,i-i ') , 


( 12 ) 


with g(j, k, l ) the Green’s function at the mesh points r m (j, k, l). This can be done efficiently 
by a three-dimensional fast convolution. The components of the electric field are determined 
as difference quotients (f.i. E c = [0 m (j, k, k, l+l)]/h z ). These components, that are 

allocated on shifted meshes (as r m (j, k, / + 1/2)), are interpolated to the Lorentz transformed 
particle positions Tj >c to find the fields E, :C in the rest frame. 

For an efficient computation of the Green’s function on the mesh, one utilizes the asymp¬ 
totic behavior for points far from the origin and calculates only near points by Eq. (J9j). 
Therefore one calculates first the anti derivative H([j — 1/2 ]h x , [k — 1/2 ]h y , [l — 1/2 \h z ) in 
one octant ( j , k, l > 0) of the near-volume and then the Green’s function itself. The field 
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in other octants follows from mirror symmetry. Technically, the change to the asymptotic 
approximation is realized with a smooth switch function 


( 


0 


if x < C 


S(x) = 1 


if x > C 2 
otherwise 


(13) 




0.5 — 0.5 cos \ ^ c 2 -Ci 


as “switched” Green’s function 


G s (r) = G(r) + (G a (r) - G(r)) S(r/ ma x(h x , h y , h z )) 


(14) 


The volume of pure asymptotic behaviour follows from the condition 


r > C 2 max(h x , h y , h z ) . 


(15) 


In Fig. [2] the Green’s function is compared with the asymptotic approximation, and the 
ratios G a /G and G s /G are plotted for cubic mesh cells {h x = h y = h z ) and for mesh cells 
with large aspect ratio (h x = 10 h y = 10 h z ). The relative deviation of the asymptotic 
function is below 0.001 for x/max(h x ,h y ,h z ) > 10. Note that the case of infinite aspect 
ratio (max(h y , h z )/h x — y 0) is very similar to the case h x = 10 h y = 10 h z in the volume 
r > 2 h x . The switch function is calculated for the parameters C 2 = 2Ci = 10. The relative 
deviation of the “switched” function is below 0.001 for the hole volume. 

C. Periodic Source and Periodic Green’s Function 

Each particle of the set S p is infinitely repeated with the spatially periodic shift r p , see 
Fig. [TJo-c. The continualization of the particles of the set gives the smooth density function 
p C)P (r) and the periodic extension is 


OO 



(16) 


The particles in S p and the continuous function p c ,p( r ) generate the periodic behavior, but 
they are not restricted to a certain volume, as the volume of one period. 

As for two dimensional charge distributions with infinite size into the third dimension, 
the potential of periodic distributions is usually not limited to a finite range. Therefore the 
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Poisson integral Eq. (J3]) diverges for all observation points of interest! It is possible to avoid 
this and to achieve convergence, by normalizing the potential 


Mr) = 0 c (r) - 0 c (r o ) = ^—J dV'x p{ r') ( 17 ) 

to be zero at a certain point r 0 . The gradient of (p c is identical to that of (j) c . In the following 
we skip the tilde and write for the potential of the periodic distribution 


0 c (r) 


1 

4tt£ 0 


OO 

dV 1 x Y p p (r' + nr p ) 

n =—oo 



47T£ f 


p OO 

J dV' x Pp(r') Y 


r + ray — r' 



The summation term in the second integral can be interpreted as Green’s function of a 
periodically repeated point charge. The Green’s function of one mesh cell follows from the 
integration 


G P (r) = 


rh x /2 r h y /2 r-hz/2 

/ / / dV ' x £ 

•J —h x /2 J —hy/2 J—h z /2 


r + ray 


|r 0 + ray - r'| 


This can be expressed by the non-periodic Green’s function as 


OO 

G p ( r ) = Y G ( r + nr p) _ G(r 0 + nr p ) . (18) 

n =—oo 


To compute the potential at mesh points k , /), caused by all mesh charges, one has 
to perform the summation 


M (j, k, l ) 


1 1 

47re 0 h x h y h z 


Y A;', /').9 P (i - j', k — k',1 — l') , 


(19) 


with q c , P {j,k,l) the integrated charge of density p CjP in mesh cell (j,k,l) and g p (j, k,l) the 
periodic Green’s function at the mesh points. Therefore the field calculation for periodic 
source distributions can be done by the same method as for a non-periodic source, only the 
Green’s function and the source distribution have to be replaced by the periodic extension 
g p and the non-repeated representation q c p . 














D. Numerical Calculation of Periodic Green’s Function 


For the numerical evaluation we set ro = 0 and split G p into the terms 

M -1 

Gi = ^2 (G(r + nr p ) - G(nr p )) , (20) 

n=—M +1 

—M oo 

G 2 = 22 (G(r + nr p ) - G(nr p )) + ^2 ( G {t + nr p ) - G(nr p )). (21) 

n=—o o n=M> 0 


The finite sum is computed directly, using Eq. (ITT)) . The number M has to be chosen 
large enough, to fulfill the condition ||r ± Mr p || > C 2 max-fh^, h y , h z } for all points in the 
computation volume, compare Eq. [15] In order to calculate the infinite sums, we use the 
asymptotic Eq. Eli 


G 2 (r) & h x h y h z 


-M 

£ 


nr r 


n\r D ) \ r + nr D 

i p/ n=M p 


h*X kyhz 


f 





? 



( 22 ) 


and the auxiliary function 


f(p,q,M > 0) 


f(p,q, 0) 


£ 

n=M 


+ 


Vv 


^(n+p) 2 + q 2 ^/(-n+p) 2 + q 2 

+ f(p,q ,!) • 


n 


|2 _i_ q2 


(23) 

(24) 


The auxiliary function f(p , q, 0) is periodic in p and can be written as fourier series 


OO 

f(p, q, 0) = flip, q) = const - 2 In \q\ + 4 £ K 0 (2nn\q\) cos(27rnp) . (25) 

n— 1 


The harmonic sum converges well for large values of |g|, due to the exp(— x)j\fx behavior 
of the modified Bessel function for large arguments. Function values for positive M are 
calculated by adding the finite sum f(p, q, M ) — f(p, q, 0). 

For small values of q the auxiliary function is calculated with help of a Taylor expansion 


f(p , q,M > 0) = f t (p , q,M) + 0 ((;p 2 + q 2 f M 10 ) (26) 
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and the Taylor polynomial 


ft(p, q, M) = (2 p 2 - g 2 ) S 3M + ^2 p 6 - 6p 2 q 2 + S 5M + 

+ ^2/ - l5pV + ^.p 2 g 4 - jjg 6 ^ S 7 . M + 

+ ( 2 / - 28A 2 + ^ pV - f pV + |^ 8 ) S 9 ,M , (27) 

with S'fc+^M = —^(k, M) /(fc!) and ip the polygamma function. To reach the desired accuracy 
it can be necessary to split f(p, q, M ) into a finite sum f(p, q, M) — f(p, q, Mf ) and an infinite 
part f(p, q, Mf) that is calculated by the Taylor approximation for a sufficiently high number 
M t . 

The numerically calculated function f(p, q, 0) and the deviation A f(p, q, 0) = | f(p, q, 0) — 
/h p (p, g, 0)| from a “high precession” result are plotted in Fig. [3] for p £ [—0.5, 0.5] and q G 
[0,1.2], The “high precession” result / hp (p, g, 0) has been calculated to numerical accuracy 
by using a Fourier series to sufficiently high harmonics or by a sufficiently high M t for the 
Taylor expansion. The highest values of A/ appear at q ~ 0.5 where the type of expansion 
is switched. The Fourier expansion is truncated after 8 terms, the Taylor expansion uses 

Mf = 16 . 


E. Mesh-Periodicity and Symmetry 

It is obvious from Fig. [T] that the choice of S p is not unique. Nevertheless the convolution 
of the continualizated density with the periodic Green’s function might result in the same 
values (pm(j, k, l ) of the potential on the mesh. This happens if the mesh supports the 
periodicity, which means the vector r p coincides with a mesh vector jh x e±i + kh y e±2 + lh z e c 

as in Fig. 0)o. 

For practical simulations it is desirable to use a mesh that supports periodicity, and to 
shift the particles of the original set by a multiple of r p so that they fill a one-period-slice 
as in Fig. |T}o and c. It stands to reason, to align the slice with respect to the mesh as in 
Fig. db. A small volume of filled mesh cells is preferable. 

Meshes that do not support periodicity can be favorable to avoid extreme aspect ratios 
of mesh cells, as it would happen for supported periodicity if the ratio of non-zero compo¬ 
nents of r p is extreme. Such a mesh, as in Fig. [4ji, gives reasonable results supposed the 
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generating set S p is sufficiently smooth. (Sudden changes of density have to be resolved by 
the discretization.) The potential has to be computed in a volume that overlaps partially 
with the shifted volume. In overlapping regions the potential is not quite identical, due to 
non-periodic discretization errors. Therefore an appropriate interpolation between poten¬ 
tials in overlapping regions (while entering one region and leaving the other) helps to avoid 
discontinuous potentials and infinite fields. 

The periodical Green’s function has point symmetry to the origin G p (j,k,l) = 
G p (—j,—k,—l). It is mirror symmetric G p (j,k,l) = G p (—j,k,l) if r p • ej_i = 0, or 
G p (j, k, l ) = G p (j, -k, l ) if r p • e ±2 = 0, or G p (j, k, l ) = G p (j, k, - 1 ) if r p • e c = 0. 


IV. APPLICATIONS 


A. FLASH Seeding Section, with 800 nm Energy Modulation of a 1.5 kA Bunch 


8] an electron bunch is energy 


In the seeding section of the FLASH free electron laser 
modulated by an external laser and send through three chicanes and several undulators 
as it can be seen in Fig. [5j To investigate microbunching effects the radiator undulator 
(after chicane 1) and the sFLASH undulators (after chicane 2) have been switched off. The 
simulation starts directly after the modulator (longitudinal position 162 m) and ends 10 m 
behind after the transverse deflecting structure TDS (at longitudinal position 204.2 m). 

The initial beam (after the modulator) has the following properties: energy £q = 585 
MeV, charge 0.3 nC, normalized emittance 1.5 /irn, peak current 1.5 kA, uncorrelated energy 
spread as = 150 keV and Gaussian longitudinal and transverse profiles. Therefore, the rms 
length is 80 fs. The beam has been energy modulated in the modulator undulator by a laser 
with wavelength A = 800 nm. The duration of the laser pulse of minimal 30 fs (FWHM) 
is short compared to bunch length, but long enough to permit pseudo periodic behavior. 
As the optical beam in the undulator is much wider than the particle beam, the energy 
modulation amplitude £ (of 250 keV) is nearly transverse offset independent. 

Effects in drift-quadrupole-lattices are the modification of particle energies, the modifi¬ 
cation of transverse optics and plasma oscillations. For the given parameters, the energy 
effect is strongest, transverse defocussing is weak and the length for a full plasma oscil¬ 
lation, in longitudinal phase space, is about L p = 120 nr [5]. Therefore the longitudinal 
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phase space distribution is nearly not altered on the short distance to the first chicane. The 
longitudinal dispersion of the first chicane R^q = 220/um couples relative energy deviations 
5s = (£ — £ 0 )/£ 0 t° longitudinal displacements 5z = RsqSs, and creates density modulation 
due to the initial energy modulation. The effect of the second chicane is weak due to its 
small dispersion R 56 = 3/im. More important is the interchange between density and energy 
modulation, as it is typical for plasma oscillations, on the distance between chicane 1 and 
3 of about Lp/6. The bunch current with its modulation is shown in Fig |6]a. Up to this 
position, the energy modulation has grown by more than a factor of ten to an amplitude of 
about 4 MeV. This amplification is beyond the linear regime and this will be exacerbated in 
the last chicane with i ? 56 = 170/irn. This longitudinal dispersion is sufficient to shear parti¬ 
cles, that were originally (during modulation) in one period, over several periods as it can 
be seen within the longitudinal phase space in Fig. [TJr. At rollover “points” in longitudinal 
phase space, where the band with particles is folded back, current spikes appear as can be 
seen in Fig. [6b. These spikes cause strong space charge fields and lead to a complicated 
energy pattern after a further drift of about 15 m, see Figs. [7b. 

Fig. [7] has been calculated by the established space charge code Astra jll] with a three 
dimensional solver based on the particle mesh method. For comparison we did non-periodic 
and periodic simulations with our implementation, called QField. The results of full-bunch 
simulations with 1E6 macro particles with Astra and QField-none-periodic can be seen in 
Figs. [T] and El showing the longitudinal phase space after chicane 3 and 15 m downstream. 
The results are in good agreement. 

To verify the periodic approach, we did a high resolution full-bunch simulation with 
20E6 macro particles with QField-none-periodic, and a periodic simulation with same macro 
particle density in only one period. The full-bunch simulations was done with longer laser 
modulation (60 fs, FWHM) to increase the length of the pseudo periodic range. Figs. [9] and 
HO] show a short part of the longitudinal phase space around the middle of the bunch, directly 
after the last chicane and 15 m downstream. Particles in gray are periodic repetitions of the 
simulated particles that are plotted in in blue. Not all simulated particles in the observation 
range are plotted to avoid large over-painted areas. (The number of plotted particles is 
equivalent for the different plots.) In the sub-panels [9b, flOh for the full-bunch, the non¬ 
periodic behavior due to finite modulation length is clearly observable. The sub-panels [9b, 
ITUb for the full-bunch simulation with longer laser modulation show a nice quasi-periodic 
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pattern that is in good agreement with the periodic results in Figs. [9k, [TOk . 


B. FLASH Seeding Section, with 266 nm Energy Modulation of a 45 A Bunch 


This example compares the three-dimensional periodic model with a one-dimensional 
impedance model, for a periodic distribution that is modulated and compressed to achieve 
current spikes of maximal amplitude. Significant self effects happen on the drift-quadrupole- 
lattice after the four magnet compressor chicane. The one-dimensional impedance model 
averages the longitudinal held of a round Gaussian beam versus the transverse offset, as 
derived in 12] . see also {3] • A periodic version of the one dimensional model was used 


for the comparison. (The periodic repetition can be realized either by a summation of the 
spatially shifted held, or by a summation of the harmonics of the impedance in frequency 
domain.) 

The setup is the same as before, but with different parameters and only one active chicane: 
beam energy £ 0 — 700 MeV, bunch charge 0.3 nC, normalized emittance 1.5 //in, peak current 
45 A, uncorrelated energy spread ag = 3keV and Gaussian prohles. The particle beam is 
energy modulated by a laser of wavelength A = 267 nm with the modulation amplitude 
£ = 0.2 MeV. The dispersion R 56 145/im of the hrst chicane is adjusted for maximum 

bunching and maximum current spikes as it can be seen Fig. [TTJ The second and third 
chicane are switched oh. 

Directly before the hrst chicane, the phase space distribution is almost not changed by self 
effects, due to the small beam current, the high particle energy and the short travelled dis¬ 
tance to this position. In contrast to the previous example, the particles are tracked through 
the magnetic helds of the chicane, but the effect of one- respectively three-dimensional self 
helds are so small, that Fig.[TT]looks identical for both cases, and it is not to distinguish from 
the corresponding hgure for a calculation without self-fields. The sharpness of the current 
spikes is not only determined by imcorrelated energy spread and longitudinal dispersion, but 
also by transverse emittance and second order coupling to the longitudinal direction. The 
current spikes are very short and the bunching factor b n = |V) exp {in2'Kz v / A) | /N is above 
ten percent even for spectral lines beyond the 30th harmonic. The required resolution of 
the longitudinal mesh has to be much better than than A/30. It is set to 0.6 nm, which 
is 450 meshlines per period of the fundamental modulation and thus below the shortest 
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possible rms length of R^cts/Sq that would appear if the particles are ideally (linear chirp, 
no non-linear effects) compressed to minimal length. There are about 220000 electrons per 
period which are directly simulated. A macro particle number equal (or higher) than the 
number of electrons can be easily handled with the periodical approach. 

The particle distribution with sharp density spikes is tracked about 23 m through 
quadrupoles and drifts to the exit of the (inactive) third chicane. Figs. \V2\ [J3] show the longi¬ 
tudinal phase space, the bunch current and the bunching for the one- and three-dimensional 
held model. The current spikes are diverged and the high harmonic content is significantly 
reduced, in good agreement of both models. It can be seen that the slope of the initial saw¬ 
tooth-like modulation has been reduced, similar to the reduction of the energy modulation 
during a plasma oscillation, after a length that is short compared to the period length. In the 
one-dimensional model all particles in one slice feel the same longitudinal held. Therefore 
the only mechanism to increase the slice energy spread is individual longitudinal particle 
motion. The longitudinal held in the spikes is strongest and causes an energy modulation 
comparable to the initial laser modulation and a large energy spread. But the space charge 
distribution significantly violates the requirements for the one dimensional approach of a 
long bunch: 707 cr±. The length of the current spike of few nano-meters times 7 is not 
large compared to the typical transverse dimension of about 100 /iin (= emittance x/3, see 
Fig. [5] with the beta function). The vectorial held of the three dimensional model depends 
on the transverse offset and is integrated along betatron trajectories. (For the other extreme 
of hat bunches 707 -C a± the longitudinal held is proportional to the transverse density, as 
dE z /dz ~ divE = p/sq.) This increases the slice energy spread as it can be seen in Fig. [13] 


C. Parasitic Heating after LCLS Laser Heater 


At LCLS a laser heater is used to increase t 
beam to counteract a microbunch instability 


re uncorrelated energy spread of the electron 


13| . Therefore the electron beam is energy 


modulated in an undulator that is in the middle of a short chicane, see Fig. [14j The periodic 
modulation directly after the undulator is (in principle) visible in longitudinal phase space. 
This is different after the dogleg comprising of the last two chicane magnets: the angular 
divergence before the dogleg, due to emittance, causes spatial longitudinal divergence due to 
the dispersion parameter r 52 . This ^-divergence is large compared to modulation wavelength 
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so that the periodic behavior is obscured, see (AS, z)-diagram in Fig. [T7j Similarly, in spatial 
space appears no pattern of the harmonic excitation, see (x, z)-diagram in Fig. |T7l Therefore 
it was unexpected to measure energy spectra with rms widths that are not proportional to 
amplitudes of the exciting laser held. The position of the spectrometer can be seen in 
Fig. [TTJ This effect is a clear indication for self effects and has been discussed, explained 
and estimated in [l(| . 

We give a simplified explanation for this effect. The transport equation from the exit of 
the laser heater undulator to an arbitrary position after the dogleg is 


x 

z 


r 11 0 
0 1 


x 0 

zo 


r 12 




5 cos (kz 0 ) 


(28) 


with x 0 , x' 0 , z 0 , 5 0 = 8cos(kz 0 ) the phase space coordinates after the undulator and x,z the 
spatial coordinates (at the arbitrary position) after the dogleg and the relative modulation 
amplitude 8. The energy spread before the undulator is neglected and for the beginning 
we neglect also the contribution of xq. Fig. [15] illustrates the transformation by Eq. [28] 
an equidistant orthogonal grid in the (x' 0 , Zo)-plane is transformed to the (x, z) plane. The 
chosen numbers are typical for the given setup. It is obvious that the original area has been 
sheared and modulated and spatial density oscillations have been created. Also particles 
with the same initial coordinate z 0 but different slope x' 0 are distributed horizontally, they 
are still located on lines. To achieve a homogeneous spatial mixing, a horizontal spread 
of xo is necessary that is coupled by rn into the horizontal direction again. Supposed the 
fluctuation of T\\Xq is sufficiently large (compared to the horizontal distance between lines 
with Ax 0 = A), then the density is uniform, as for instance in Fig. [171 This is different in the 
vicinity of zero crossings of r xl as in Fig. [16] with r n and optical functions (3 x ,/3 y . At these 
positions appear spatial density fluctuations and space charge fields cause an additional 
energy modulation, which is correlated with the original modulation. 

Figs. [TT] [18] and [T9] are results of a numerical simulation with the periodic space charge 
solver of the setup described in [10] (S = 135 MeV, <j£ = 2 keV, peak current 37 A, charge 
250 pC, normalized emittance 0.4 /mi, modulation wavelength A = 758 nm). The beam is 
heated to a £ = 8 keV. The blue particles in the (AS, z)- and (x, z)-diagrarris are particles 
that have been in the same longitudinal period during modulation, the gray particles are 
a periodic repetition. Fig. [17] refers to a position short after the laser heater chicane. The 
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particles behave as expected by a theory neglecting self effects: the oscillating behavior is not 


visible, the particles are longitudinally smeared over a length larger t 


ran the period A. The 


spectral density agrees with the predicted double bump structure Q. Fig. ITT! refers to the 
position 15.5 m, close to the zero crossing of rn. In the spatial diagram the periodic pattern 
appears as microbunching, although the relative longitudinal position of all the particles 
is not changed. The space charge fields due to the micro structure changes the particles 
energy and alters the energy distribution: the spectrum is different and its rrns width is 
increased. At the position 17.5 m shown in Fig. [T9J no microbunching is obvious but the 
energy distribution is significantly widened to ~ 12 keV. Fig. I2U1 shows the rms energy 
spread along the beamline for different heating strengths. The blue 2 keV curve is calculated 
without heating, the thick black line refers to a heating (from initially 2 keV) to 8 keV as in 
the previous figures. The parasitic heating effect can be clearly seen. The characteristic of 
this heating described in Fig. [21] showing the rms energy spread at the spectrometer entrance 
versus the initial energy spread. It is remarkable that parasitic heating changes not only the 
rms width but also the shape of the spectrum, see sub-panels with spectra in Figs. [171 [T8] 
and m 

This calculation has been performed with 1E6 macro particles per period although the 
number of electrons per period is 0.55E6. The longitudinal resolution is Az = A/50, the 
transverse resolution is Ax = A y = 'yAz, so that the Lorentz transformed mesh has the 
same resolution in all directions. The phase space pictures have been computed separately, 
with a resolution of 20 meshlines per period. 

Fig. [22] shows the published measurement 10] together with a periodic and non-periodic 
simulation with QField. Although the shapes agree qualitatively, the results are not in good 
agreement: (a) the measured heating is stronger (up to 28 keV where 10 keV are expected 
without additional heating), (b) the energy values without increased heating are different 
(here: ~ 11 keV, 20 keV, measured: ~ 16 keV, 28 keV). 


V. CONCLUSION 

In this paper we developed a method to calculate the potential and electro-magnetic held 
of a periodic charge distribution that moves in free space with constant velocity. The charge 
distribution is three dimensional and periodic in one direction of space. The numerical 
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technique is a modification of the well known particle-mesh-method with fast convolution of 
the discretized source with the Green’s function of one cell. An efficient method has been 
described to compute the periodic Green’s function so that only particles, representing one 
period, have to be considered for the convolution. 

The method has been demonstrated with examples of pure spatial periodicity into the 
direction of motion. For example (A) a comparison with non-periodic simulations is shown. 
Although the quasi-periodic part of the bunch is quite short, a good agreement is achieved. A 
non-periodic simulation for example (B) was beyond the available computational capacities. 
For example (C) some points of the curve in Fig. \T2\ > have been verified with large effort by 
non-periodic simulations with lower spatial resolution. The CPU time requirement for each 
point was about 1000 times larger despite the lower resolution. 

Not all capabilities of the method have been demonstrated by these examples. For in¬ 
stance a density modulated beam with energy chirp represents a periodicity in geometric- 
and momentum-space. Due to longitudinal and transverse dispersion the periodicity vector 
r p changes its length and direction and deviates from the direction of motion. This happens 
in bunch compressor chicanes. 

The periodic approach is numerically much more efficient than non-periodic simulations 
of bunches with quasi-periodicity, as the number of particles and the volume of the problem 
domain is reduced by orders of magnitude. 
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FIG. 1: Periodic particle distributions and periodic shift vector r p : (a) pseudo-periodic random 
distribution with periodic behavior; (b) generating set S p and periodically repeated particles in blue 
resp. gray; (c) equivalent generating set in slice volume between two planes shifted by r p ; (d) other 
equivalent set in different volume slice; (e) equivalent generating set without sharp truncation. 
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FIG. 2: The Green’s function G, the asymptotic approximation G a and the “switched” approxi¬ 
mation G s for switch parameters C 2 = 2Ci = 10: (a) and (c) comparison of G and its asymptotic 
approximation G p ; (b) and (d) ratios G a /G and G s /G; (a) and (b) for cubic mesh cells; (c) and (d) 
for mesh cells with large aspect ratio (h x = 10 h y = 10 h z ). The Green’s funtion is only calculated 
for discrete mesh points (x = nh x ). 
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FIG. 3: Function f(p,q, 0) and absolut error Af(p,q,0) of the numerical calculation. As / is 
singular at the origin (p = q = 0), it truncated for the plot to values below 10. The transition 
from the Taylor approximation (for small \q\) to the Fourier series is at |q| = 0.5. The Taylor 
approximation is used for the upper part of the sum from Mf = 16 to infinity. The harmonic 
expansion is truncated after 8 terms. 
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FIG. 4: Equivalent generating sets S p for the same periodic distribution: (a) volume is not limited 
to a one-period-slice; (b) particles in one-period-slice and mesh that supports periodicity; (c) 
particles in equivalent one-period-slice; (d) particles without sharp truncation and mesh that does 
not support periodicity. 
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FIG. 5: Seeding section of FLASH: optical functions and layout. 



(a) (b) 

FIG. 6: Application A: current (a) before and (b) after the third chicane, for 30 fs (FWHM) laser 
modulation, calculated with 20E6 macro-particles. 
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FIG. 7: Application A: longitudinal phase space (a) direct after chicane 3 and (b) 15 m behind 
for 30 fs (FWHM) laser modulation, calculated with 1E6 macro-particles, by Astra. 



FIG. 8: Application A: longitudinal phase space (a) direct after chicane 3 and (b) 15 m behind 
for 30 fs (FWHM) laser modulation, calculated with 1E6 macro-particles, by QField. 
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FIG. 9: Application A: longitudinal phase space direct after chicane 3. (a) Calculated with QField- 
non-periodic, 1E6 macro-particles, 30 fs (FWHM) laser modulation, (b) calculated with QField- 
non-periodic, 20E6 particles, 60 fs (FWHM) laser modulation, (c) calculated with QField-periodic. 
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FIG. 10: Application A: longitudinal phase space 15 m behind chicane 3. (a) Calculated with 
QField-non-periodic, 1E6 macro-particles, 30 fs (FWHM) laser modulation, (b) calculated with 
QField-non-periodic, 20E6 particles, 60 fs (FWHM) laser modulation, (c) calculated with QField- 
periodic. 
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FIG. 11: Application B: longitudinal phase space, current and bunching after 1st chicane. 




FIG. 12: Application B: longitudinal phase space, current and bunching after 3rd chicane, calcu¬ 
lated with one dimensional impedance. Chicanes 2 and 3 are off. 




FIG. 13: Application B: longitudinal phase space, current and bunching after 3rd chicane, calcu¬ 
lated with three dimensional impedance. Chicanes 2 and 3 are off. 
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FIG. 14: Application C: the LCLS injector layout showing laser heater, transverse RF deflector 
OTR/YAG screens wire scanners, and spectrometers, from [10] 
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FIG. 15: Application C: visualization of transport transformation of (zq,x'q), after laser heater 
undulator to (x,z), at 15.5 m position for particles without offset xq. 
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FIG. 16: Application C: beta functions and rn from the end of the laser heater undulator to 
an arbitrary lattice. The plot range is from the exit of the undulator to the beginning of the 
spectrometer magnet. 
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FIG. 17: Application C: longitudinal phase space (A £,z), top view (x,z) and energy spectrum for 


the position 11 m, short after the chicane. 
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FIG. 18: Application C: longitudinal phase space (A S,z), top view (x,z) and energy spectrum for 
the position 15.5 m, near the zero crossing of r\\. 
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FIG. 19: Application C: longitudinal phase space (A £,z), top view (x,z) and energy spectrum for 
the position 17.5 m. 
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FIG. 20: Application C: energy spread after the laser heater chicane as function of the beamline 
coordinate. The initial spread varies due to different laser modulation. 



FIG. 21: Application C: energy spread at the spectrometer (position 19 m) as function of the 
spread directly after the laser heater chicane (position 11 m). The initial spread varies due to 
different laser modulation. 
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FIG. 22: (a) Measurement 



(b) simulation with QField-periodic and QField-non-periodic. 
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